arXiv:1505.06937v2 [physics.soc-ph] 15 Sep 2015 


Phase transitions in Pareto optimal complex networks 


Luis F. Seoane^’^ and Ricard Sole^’^’^ 

^ ICREA-Complex Systems Lab, Universitat Pompeu Fabra - PRBB, Dr. Aiguader 88, 08003 Barcelona, Spain 
^ Institut de Biologia Evolutiva, UPF-CSIC, Passg Barceloneta, 08003 Barcelona 
® Santa Fe Institute, 1399 Hyde Park Road, New Mexico 87501, USA 

The organization of interactions in complex systems can be described by networks connecting 
different units. These graphs are useful representations of the local and global complexity of the 
underlying systems. The origin of their topological structure can be diverse, resulting from different 
mechanisms including multiplicative processes and optimization. In spatial networks or in graphs 
where cost constraints are at work, as it occurs in a plethora of situations from power grids to the 
wiring of neurons in the brain, optimization plays an important part in shaping their organization. In 
this paper we study network designs resulting from a Pareto optimization process, where different 
simultaneous constraints are the targets of selection. We analyze three variations on a problem 
finding phase transitions of different kinds. Distinct phases are associated to different arrangements 
of the connections; but the need of drastic topological changes does not determine the presence, 
nor the nature of the phase transitions encountered. Instead, the functions under optimization do 
play a determinant role. This reinforces the view that phase transitions do not arise from intrinsic 
properties of a system alone, but from the interplay of that system with its external constraints. 

PACS numbers: 64.60.aq, 64.60.Bd, 64.60.-i, 87.55.de 


I. INTRODUCTION 

Optimization is a key goal in engineered systems and is 
traditionally assumed to be part of the intrinsic dynamics 
of natural evolving systems [T] . The engineering perspec¬ 
tive, associated to man-made objects and structures, is 
specially obvious when dealing with large-scale, intercon¬ 
nected units, as it occurs in very large integrated circuit 
design Hi] or spatially-extended infrastructures such as 
power grids [S? -[7] and transportation or distribution 
networks jSHia. In these cases, engineers cope with in¬ 
terfering constraints related to materials, space, packing, 
wiring, or dissipation costs. The staggering complexity 
of these designed systems can be addressed by algorithms 
that deal with multidimensional problems. 

In biological systems, important network topologies 
have been shown to result from optimality [131 mi . These 
include transportation networks in living organisms Ha¬ 
ul] where optimization is reached by means of fractal 
trees that guarantee a low cost and efficient location of 
resources. Similarly, neural circuits display optimal fea¬ 
tures over a wide range of scales H [20U24] . The packing 
and interconnectivity in some cortical areas seem com¬ 
patible with design principles shared by high-density elec¬ 
tronic designs |3]. 

In all the previous examples tradeoffs between effi¬ 
ciency and cost are present. Packing many components 
in a given spatial domain is desirable because of cost 
minimization of connections, but dissipation of energy or 
wiring constraints will also be at work. What kinds of 
topologies result when considering multiple constraints? 
This problem has been addressed by explicitly introduc¬ 
ing efficiency measures E (such as average path length) 
along with cost constraints C (such as number of connec¬ 
tions of a given graph) [25U27] . A similar example in an¬ 
other field models languages as a network of associations 


between objects and words, and considers language evo¬ 
lution through a least effort process PHHSO] . Here, the 
cost-efficiency conffict is mapped onto coding/decoding 
efforts for users of an economic (while ambiguous) lan¬ 
guage. Both in the linguistic and the network examples, 
the tension between opposite demands leads to phase 
transition phenomena. We can wonder, in a more gen¬ 
eral note, when and how will phase transitions arise dis¬ 
regarding of the details of the problem in hand. 

Within the context of network optimization, we con¬ 
sider the set r of all connected networks 7 S T involving 
N nodes and any number of links. Latter on we will de¬ 
fine an efficiency {E{j)) and a cost {C{j)) based on the 
structure of each network 7, so that a series of optimiza¬ 
tion problems can be posed. For each such graph, we can 
also introduce a global energy function ^(7) that takes 
into account our optimization goals. The most straight¬ 
forward way to do this is through a linear combination: 

L!(7,A) = Aii;(7) + (l-A)C'(7), (1) 

with A S [0,1] a tunable parameter that weights the im¬ 
pact of each contribution. This is precisely the strategy 
in previous accounts of these problems [25H3T] . If A = 1 
only efficiency constraints will be at work, whereas A = 0 
would ignore this component. 

Such a global energy function results in very illustra¬ 
tive visualizations of the optimization process through 
the notion of a potential landscape. Assuming minimiza¬ 
tion, more optimal network architectures lay deeper in a 
potential well when we plot ^(7, A) for every 7 superim¬ 
posed on an arbitrary network morphospace, as in Fig. 
Sa). Note, however, that this is a limited picture: a fixed 
value of A is necessary to generate one potential land¬ 
scape. Changing the parameter modifies the landscape 
rendered by ^(7, A) and, accordingly, the underlaying op¬ 
timization problem. To achieve a more general compre- 


2 


hension we should not only allow scenarios with different 
values of A, but we must also question the hypothesis of 
linearity introduced by equation §• Therefore, we con¬ 
sider Pareto (or Multi Objective) Optimization 
whose solution is not a global optimizer dependent on 
external parameters (i.e. the absolute minimum of some 
potential landscape), but a collection of solutions that 
attempt to satisfy both (cost and efficiency) constraints 
simultaneously. These Pareto optimal designs constitute 
the so-called Pareto front: the most optimal tradeoff pos¬ 
sible between the targets involved (Fig. Be)). 

In this paper we look at such tradeoffs for a series of 
optimization problems defined upon complex networks. 
How the different constraints are satisfied depending on 
a series of factors is interesting in itself; but if we re¬ 
late the Pareto optimal designs to the linear optimiza¬ 
tion problem posed by equation Q , a richer phenomenol¬ 
ogy unfolds. This stems from a deep connection between 
the Pareto front and phase transitions and other impor¬ 
tant features of thermodynamic theory |35] . In this con¬ 
text, the energy landscapes that we can generate through 
equation Q often result very useful to complete the pic¬ 
ture. 

The paper is organized as follows. Sec. [H] provides the 
basic formalism of Pareto optimization and the connec¬ 
tion between Pareto optimality and statistical mechanics 
[35j . This theory homes in phase transitions for the prob¬ 
lems investigated and is the framework used to analyze 
most of the results obtained. Also in Sec. we pose 
a series of Multi Objective Optimizations upon complex 
networks whose solutions are detailed in Sec. Ell and 
discussed, with closing remarks, in Sec. |IV| Appendix [A| 
explains some numerical aspects of the current work. 


II. PHASE TRANSITIONS IN THE PARETO 
FORMALISM 

A recent contribution proposed a multiobjective opti¬ 
mization approach to statistical mechanics that gener¬ 
alizes key concepts in thermodynamics for any Pareto 
optimality problem [35]. The fundamental connection is 
between the Pareto front and the Gibbs surface. This 
surface (defined through the thermodynamic potential 
G = G{U,S,V)) is linked to the equilibrium state of a 
thermodynamic species. Then, its concavities and non- 
differentiable edges underly the existence of first and sec¬ 
ond order phase transitions [MllE]. 

A similar point can be made for systems that opti¬ 
mize a set Tf of target functions Tf = {ti ,..., tx}, whose 
Pareto front encodes phase transitions and critical points 
in its shape. This connection between Pareto optimal 
systems and thermodynamics is explained at length in 
[55] . There, general Pareto optimal designs (not nec¬ 
essarily networks) that minimize an arbitrary number 
AT > 2 of targets is discussed. Here we sketch the theory 
with just K = 2 optimization targets and using con¬ 
nected networks with fixed number of nodes (7 S P). 
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FIG. 1 : A two dimensional example of Pareto opti¬ 
mality. (a) 7 £ r are all possible connected networks 
with a given number of nodes. They populate some network 
morphospace where we seek those graphs minimizing some 
measurable feature. If we deal with just one fitness function, 
an energy landscape can be defined and the optima are easily 
found at the bottom of energy wells, (b) If more than one op¬ 
timization target are at play, this landscape picture falls apart 
and we need to adopt a Pareto optimization approach. Then 
our task is to find a set of Pareto optimal solutions (Hr C P) 
that minimizes all targets (here ti and t2) simultaneously. 
These functions map each network 7 G P into R^. The subset 
of Pareto optimal solutions is mapped into the Pareto front 
(thick gray curve). Along this curve it is not possible to im¬ 
prove both ti and t2 at the same time, (c) It is convenient 
to introduce Pareto dominance. The black circle represents a 
net with lower ti and t2 than (thus it Pareto dominates) those 
within the gray square. That same filled circle is dominated 
by those networks projected between the dashed lines and the 
Pareto front, (d) Recovering a single objective optimization 
through a linear combination of the targets (Eq. 0 ) is akin 
to choosing a direction in R^ and seeking the extreme of the 
front along that direction. 


This suffices to illustrate the relevant aspects of the theo¬ 
retical framework. We remit the reader to |35j for further 
details. 

Henceforth, Tf = {ti,t 2 } are any two real valued func¬ 
tions that we can measure on any network. These func¬ 
tions project each network into Tf{’j) = (^ 1 ( 7 ), t 2 ( 7 )), a 
point in the plane that we term target space (Fig. B^)- 
(b)). In that plane we can solve the Multi Objective 
Optimization (MOO) problem consisting of the simulta¬ 
neous minimization of both targets. Therefore we define 
Pareto dominance (Fig. B*^))- ^ network 7 ^; S F domi¬ 
nates another 7 ^ G F (and we note it jx ^ "fy) if 

tkilx) < tkily) V fc = 1,2; 

3fc'e{I,2} I tk'i'^x) <tk'{ly). (2) 

In this case ^x is objectively better than 7 j, and we can 
dismiss the later. It is often the case that pairs of net- 
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works are mutually non dominated. Then we cannot 
choose between them without introducing a bias towards 
either ti or t 2 - We have to avoid this bias to find Pareto 
optimal designs. We say that a network 7 ^; S P is Pareto 
optimal if it does not exist any other 7 j, G P such that 
Hence 


Hr = {7a: e P I ^ 7j; G P, 7j; ^ 7^:} (3) 

is the set of all Pareto optimal solutions, which consti¬ 
tutes the solution of the MOO problem. 

The set Hp C P is projected onto the target plane 
through Ty(nr) where it represents a limiting frontier of 
the whole T/(r), known as the Pareto front. For K = 2 
this frontier is a curve that implements a bijective func¬ 
tion of <1 and ^2 (Fig- Qc)). Hp is not a standard global 
optimizer. Instead, it comprises a collection of valid net¬ 
works that embody the optimal tradeoff between many 
targets such that, as we move through it, we cannot im¬ 
prove ti without worsening ^2 and vice-versa. 

We can now add a further demand that a global energy 
function 


11 ( 7 ; A) = Ati(7)-k (1 - A)t2(7) (4) 

be minimized. (This equation is equivalent to Eq. [^with 
the more general ^ 1 , 2 ( 7 ) instead of efficiency and cost. 
We rewrite it here for convenience.) Note that this is the 
simplest Single Objective Optimization (SOO) that we 
can built with Tf. Indeed, each possible A G [0,1] poses a 
different SOO problem, potentially with a distinct global 
solution. We check now how solutions of this SOO family 
are related. 

By setting fix a value of A we choose a direction in 
the ti — O plane along which the minimization proceeds 
(Fig. Hd)). From Eq. networks projected onto a 
straight line with slope d = —A/(l — A) have the same 
global energy Q. Networks on straight lines pushed fur¬ 
ther against the front present lower energies, until one 
solution 7 a G Hp is singled out as global minimum of 
It usually lays where the slope of the front is pre¬ 
cisely d = —A/(l — A). If the front Hp is convex and its 
derivatives are well defined, we can smoothly sample the 
front by slowly tuning A. Any appropriate measurement 
that we perform on 7 a will be a smooth, differentiable 
function of A itself. 

If a cavity exists in T/(np), the global solutions bypass 
those that lay inside (Fig. [^a)). For a certain value 
A'^, two very different solutions (y^c and 7 ac) located far 
apart in the ti — t 2 plane are simultaneously optimal. 
Solutions between 7 ^^ and 7 ^^ (inside the cavity) never 
get to be global optima. For A < A'^ we remain at one side 
of the cavity where a smooth sample of the Pareto front is 
still possible. The same happens for A > A'^ at the other 
side of the cavity. But at the characteristic value A'^ there 
is a sudden jump between and y^c and, consequently, 
in any order parameter that we can measure on 7 a when 
plotted as a function of A (Fig. ib)). This is a first 
order phase transition. A similar transition happens if 





FIG. 2 : Phase transitions in Pareto optimal systems. 

(a) A first order phase transition takes place in MOOs whose 
Pareto front presents a cavity. Solutions of the convex hull 
of the front are SOO global optima for a range of A, while 
those inside the cavity are not. (b) This renders a gap in 
any adequate, measurable property of the global solutions, 
(c) Second order phase transition are associated to ill-defined 
derivatives - i.e. sharp edges in the front. Those designs at 
the edge are persistently optimal foe a range of the control 
parameter A. (d) Hence, any property of the global solution is 
constant as a function of that same control parameter for the 
designated interval, (e) Second order phase transitions might 
also take place at the extremes of the Pareto front. We say 
then that the front has an abrupt ending (or that it ends in 
an abrupt manner) as in the lower-right extreme of this front, 
(f) Any order parameter presents the well-known functional 
dependence of second order phase transitions. 


the Pareto front is a straight line (then A'^ = —d°/(1—A^), 
where is the slope of the front) with the singularity 
that all Pareto optimal solutions are also global optima 
for A = A^. 

If the front is convex but its slope is ill-defined any¬ 
where (Fig. [^c)) there is a range of values A G (A“, A+) 
for which the global optimum remains the same. Because 
the optimum does not change within this singular range, 
if we plot any order parameter its derivative (with respect 
to A) will be zero. The same derivative will be different 
from zero anywhere else. This results in a characteristic 
plot for any order parameter (Fig. id)) associated to 
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second order phase transitions. 

Such non-differentiability of the front might be shifted 
all the way towards one of its extremes. For A G (0,1) 
we seek global optima where the slope of the front is 
any —A/(l — X) => d G (—oo,0). Well-behaved Pareto 
fronts should roll smoothly with slopes that completely 
span the d G (—oo, 0) interval; otherwise the front will be 
exhausted while more SOO problems can still be defined. 
If these slope values are not exhausted we say that the 
Pareto front ends in an abrupt manner, meaning that its 
slope as we approach its leftmost end is some d*_ > —oo, 
or its slope as we approach its rightmost end is some other 
d*^ < 0. These abrupt endings of the front will be optima 
already for some A* = —d^/ (1 — d^;)- We can still define 
SOO problems for A beyond these limits, for which the 
solutions at the corresponding extremes of the front will 
be persistently optimal (Fig. [^e)). This implies again 
a discontinuity in the derivative of any order parameter, 
the fingerprint of second order transitions (Fig. [^f)). 


A. Multiobjective optimization of complex 
networks 


Complex graphs are a great testbed to illustrate this 
theoretical framework. They allow us to define prob¬ 
lems of increasing difficulty where first and second order 
transitions arise. Seminal work on network optimiza¬ 
tion addressed the problem from an SOO perspective 
[TTl [23 Uni [311 [35], so some of our results can be put 
in context. Another advantage of complex networks is 
that good optimizers can be produced in the computer, 
simplifying the empirical work. 

We propose three problems based on the conflict be¬ 
tween the average path length between nodes and the 
density of edges, which roughly inform us about diffu¬ 
sion efficiency [33] (for which low average path length 
is desired) and implementation costs (lower for sparser 
networks). Consider first the topological (or standard) 
average path length: 


(0‘(7) 


1 


(5) 


where dX(^) denotes the distance (in number of edges) 
between nodes ni,nj G 7 along the shortest path that 
connects them; and the topological (or standard) link 
density: 


/w = (0) 

P ij 

where the adjacency matrix A( 7 ) = { 0 ^( 7 )} presents 
aij{j) = 1 if two nodes are linked in 7 and aijipf) = 
0 otherwise. and .Z* are normalization constants 

discussed below. 

The superindices in (1)* ( 7 ), ^*( 7 ) indicate that we deal 
with the topological (or standard) average path length 


and link density, in which edges cost 1 unit. Geometric 
costs can be included if nodes are distributed, e.g., over 
a Euclidean space. Let d® ( 7 ) be the Euclidean length 
of the shortest path connecting Ui and Uj in network 7 
- i.e. the sum of the Euclidean lengths of the edges in 
the shortest path between these nodes provided that 7 
is embedded in some geometric space. We introduce the 
geometric (or weighted) average path length: 


(0®(7) 


d) 


2 


(7) 


The shortest Euclidean distance possible between two 
nodes lij{^) only enters equation Q if a direct link be¬ 
tween Ui and Uj is present in 7 (in that case dfj{'j) = 
hjil))- This kji'j) allows us to introduce the geometric 
(or weighted) link density: 

P»(7) = (8) 

^ id 

Just as before, (1)® ( 7 ) and pP{j) (note the su¬ 
perindexes indicating their geometric dependence) are 
normalized by and Z^. A clique, or fully connected 
network ( 7 c), has the shortest average path length 
possible always. As we will see later, this means that 70 
is Pareto optimal always, so we base our normalization 

on it: Z*/)® = Ei 7 C^®(bi; 7 c)/ 2 , Zj, = E*j ay( 7 c)/ 2 , 
and Zl = 

We combine (1)*^® ( 7 ) and p*'^®( 7 ) as targets in different 
ways to generate three MOO problems: 

(A) Fully topological problem, with ti = (/)*( 7 ) and 

^2 = Note that the geometry does not play 

any role in this case. 

This version was originally studied in [23126| from 
an SOO perspective. From that approach, only the 
clique and star graphs appear relevant (as discussed 
in [23) as the representatives of two phases at ei¬ 
ther side of a discontinuous phase transition. We 
show how this fits parsimoniously within the frame¬ 
work presented in [33. But besides, we discuss 
now the whole Pareto front - its relevant shape and 
some of its constituents. This front includes non¬ 
trivial complexities well differentiated from the star 
and clique, and it presents connections to critical 
systems discussed elsewhere [40] . 

(B) Partly geometrical problem, with ti — (1)* ( 7 ) (the 
same as above) and ^2 = P® ( 7 ) ■ Geometry, through 
t 2 , plays a relevant role now. Since the disposition 
of the nodes in space matters, we study this MOO 
in two different cases: i) nodes scattered randomly 
over the [ 0 , 1 ] x [ 0 , 1 ] square in and ii) nodes 
spaced evenly over a circle of radius 1 . 

In this problem we still use the topological average 
path length, meaning that we seek to minimize the 
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number of hops or the number of relay stations be¬ 
tween arbitrary pairs of nodes. To think about this 
problem we can picture an infrastructure such as a 
subway network whose contractor wishes to mini¬ 
mize the length of line built while the users want 
to avoid transfers between lines. 

(C) Fully geometrical problem, with t 2 = {l)^ ( 7 ) and 
t 2 = p®( 7 ). In this case the geometrical cost is 
important for all targets involved. Again, the dis¬ 
position of the nodes matters and again we study: 
i) nodes scattered randomly over the [ 0 , 1 ] x [ 0 , 1 ] 
square in and ii) nodes spaced evenly over a 
circle of radius 1 . 


III. RESULTS 

Three relevant topologies indicate major feats of all 
our Pareto fronts. The most prominent one is the clique: 
a fully connected network that presents the largest num¬ 
ber of links possible, thus maximizes edge density and 
minimizes the average path length always. This guaran¬ 
tees that the clique is always Pareto optimal. It marks 
the top-left boundary of the Pareto front, as illustrated 
in Fig. §a). This is true for all problems considered in 
this paper. 

The star presents a hub to which all other nodes are 
connected, while non-hubs are not connected to each 
other. There are N possible star graphs. If geometry 
is not considered, all of them are equivalent. When ge¬ 
ometry intervenes and nodes are spaced over a circle all 
N stars are equivalent as well. All possible trees consist 
of as many edges as the star but, if geometry matters, 
only the minimum spanning tree (MST) minimizes al¬ 
ways the edge density (^ 2 )- The MST is Pareto optimal 
whenever geometry is relevant and it always indicates 
the end of the Pareto front opposite to the clique (at its 
bottom-right). 


A. Fully topological problem 

This case has been studied as an SOO through equa¬ 
tion Q [33 [3^. That solution is incomplete from an 
MOO perspective which was not the chosen paradigm in 
those works anyway. This problem has the advantage 
that its front (Fig. [^a)) can be found analytically and 
the phase transitions derived from it are independent of 
the number of nodes. We cannot guarantee the same for 
the variations studied later. 

Because we normalized both targets using the clique 
as a reference, this network is mapped into (1,1) in the 
ti — t 2 plane. Any graph will have less edges than a 
clique, thus the set F of all connected networks lays be¬ 
low t 2 = 1 in the target plane. The lower boundary 
of t 2 is achieved by connected networks with the min¬ 
imum amount of edges possible (TV — 1). There are 



FIG. 3 : Pareto front of the fully topological problem. 

(a) The front (solid gray curve) is a straight line connect¬ 
ing two phases: a star and a clique. The slope of the line 
(F = —1 determines that at = 1/2 a first order phase 
transition takes place. All networks laying on the front are 
global SOO optima at that critical value. Among them we 
find networks produced by attaching links to a star and oth¬ 
ers radically different from the star and from the clique (note 
the two graphs marked A and B: only one of them can be 
produced by attaching edges to the star). (b) All core graphs 
for have been listed for N < 5 . Beyond that, it becomes in¬ 
creasingly difficult to count how many there are or even to 
tell apart two different ones. 


a collection of such graphs, from the star to a linear 
chain - in between lay all possible trees. All them have 
= 1 /AT, which tends to 0 as A goes to infinity. 
The average path length of these networks varies between 
that of the star (2(iV — 1)/TV —)■ 2) and the linear chain 
((AT^ — 1)/3(AT— 1) —>• -boo). These minimally connected 
graphs lay on a horizontal stretch of the ti — t 2 plane 
(dashed line in Fig. ia)). 

Among these trees (all with the same ^2 = 1/AT), the 
star is the one with the lowest average path length, hence 
it is Pareto optimal. Any other network with a lower ti 
must have more links than the star, the clique setting 
the lower t 2 bound. Thus the Pareto front must lay on a 
curve connecting the clique and the star - i.e. connecting 
( 1 , 1 ) and ( 2 , 0 ) in the U — t 2 plane. 

We appreciate the following facts: i) The edge density 
is a function of the number of links alone and it does 
not depend on the topology of the network, ii) Given a 
network that is Pareto optimal, we generate new Pareto 
optimal networks by simply adding new connections. As 
an instance, the star is Pareto optimal and all its nodes 
are 1 edge apart from the hub and 2 edges apart from 
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FIG. 4: Partly geometrical problem on nodes scattered over a plane. (a) The front follows the archetype of the 
topological problem with two roughly perpendicular stretches that trade off between the clique (top-left), the star, and the MST 
(bottom-right). Incomplete cliques are reached after a second order phase transition because the Pareto front ends abruptly in 
its top-left (inset). The other extreme of the front ends smoothly, (b) A sharp edge indicates a second order phase transition 
with the star graph being optimal for a range A € (A“,A^). (c) Plotting an order parameter as a function of A reveals both 
transitions at A ~ 0.61 and at A“ ~ 0.01 and A^ ~ 0.3 (inset). The star is optimal in the range A G (A“, A^), thus any order 
parameter is constant in that range, (d) The non-analyticity of the Pareto front is inherited by the energetic landscape also 
as a sharp edge. The SOO is vividly illustrated thanks to this potential landscape, whose minimum is occupied by one same 
network for several values of A. As lambda changes, the potential well is deformed until the minimum drifts away from the 
sharp edge. 


each other. Then, new edges can only be added that 
connect directly two non-hub nodes, turning a distance 
d^ihj) = 2 into d}{i,j) = 1; but not affecting the network 
in any other respect. Put otherwise, once a network is 
Pareto optimal any addition of links has got only local 
effects in its average path length. 

Adding new links to the star results in more Pareto op¬ 
timal networks, the number of which grows combinato- 
rially (that scaling saturates as we approach the clique). 
Take apart the iV — 1 non-hub nodes of a star: any net¬ 
work that we implement on this subset of nodes (con¬ 
nected or not), and which is subsequently embedded on 
the original star graph through the hub, is Pareto opti¬ 
mal. It is a sujficient (but not necessary) condition for 


a network to be Pareto optimal to contain a hub (Fig. 
I^a)). The necessary condition for Pareto optimality is 
that every node is at maximum 2 edges apart from each 
other. 

From any Pareto optimal network (with or without a 
hub), adding new edges always generates new Pareto op¬ 
timal graphs. Repeating this operation we always reach a 
clique, but this process is not reversible: Take the clique 
and delete connections randomly with the condition that 
your network remains Pareto optimal after every dele¬ 
tion. No rearrangement of the edges is allowed. Let this 
process continue until we cannot remove any link with¬ 
out violating the Pareto-optimality condition. This algo¬ 
rithm might yield a star or any other graph from a col- 
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lection of irreducible Pareto optimal networks, which we 
call Pareto core graphs. The star is the core graph with 
less edges possible. We can only construct these networks 
as described, since other defining regularities are not ap¬ 
parent - beyond the optimality condition that every node 
is at most 2 hops away from each other. Some of these 
graphs are represented in Fig. ib) for N = 2,... ,7. 
The complexity scales from 1 core graph for N = 2,3; 
to two core graphs for TV = 4; to three for TV = 5; to 
an unknown number for TV > 6. For larger TV it also 
becomes increasingly difficult to determine whether two 
core graphs are the same, given their invariance under the 
labeling of the nodes. Note that core graphs are Pareto 
optimal. They are representative of the staggering com¬ 
plexity contained in the front (which grows combinato- 
rially) and they cannot be trivially composed as a mix¬ 
ture of stars and cliques. Because of this they constitute 
Pareto and global optima that have not been previously 
reported. 

Even if we cannot list down all Pareto optimal net¬ 
works, we can find where they live on the ti — t 2 plane. 
Adding one edge always modifies (/)* ( 7 ) by an amount 
A (T)* = —1/TV(TV — 1), thus ti is decreased. The same 
operation increases ^2 by Ap‘ = 1/TV(TV — 1). Because 
Ap*/A{l) = —1 does not depend on the number of 
edges, Pareto optimal graphs thus generated lay on a 
straight line with slope = Ap*/A{l) (Fig. ia)). 
Such a front implies a first order phase transition at 
A'^ = —d‘^/(l — d'^) = 1/2. The clique and the star are 
found at either phase (correspondingly for A > A'^ and 
A < A*^). Right at the critical value A"^ any Pareto opti¬ 
mal network is a global optimum. The degenerated com¬ 
plexity at A = A'^, with so many and structurally different 
optimal solutions, presents interesting connections with 
critical phenomena and neutral theory that we explore 
in future work |40j . The plot of any order parameter as 
a function of A (not shown) just presents a gap between 
two constant values. 


B. Partly geometrical problem 

Figure l^a) provides an archetype for the Pareto front 
that will be repeated (with variations) in the more elab¬ 
orated MOO problems. Our fronts will present a first, 
stepped stretch that trades off between the clique (top- 
left) and some intermediate networks (usually the star); 
and a second, flat stretch with little variation in the ver¬ 
tical dimension (^ 2 ) and a broad variation in the hori¬ 
zontal axis (ti). In the previous case, this second stretch 
(dashed line in Fig. [^a)) does not belong to the front, 
but it will in the following problems. 

1. Nodes scattered over a plane 

Figure]^ a) shows the first example of this archetype. 
A very stepped stretch of the front trades off between 


the clique and the star just as before. However, this is 
a convex curve now, which we discuss below. The sec¬ 
ond archetypal stretch of the front trades off between the 
star and the MST, and is mapped onto an almost hori¬ 
zontal curve in the ti — t 2 plane with a slight convexity. 
This Pareto front ends up smoothly in its bottom-right 
extreme, so we dismiss any phenomenon associated to 
it. Because the whole front is convex first order phase 
transitions are ruled out. 

The first, stepped stretch of the front (Fig. |^a), inset) 
presents a feature that appears in most subsequent cases. 
This stretch is a convex curve that ends abruptly (with 
the notion of abruptness introduced before). This indi¬ 
cates that a second order phase transition takes place. 
This transition trades off between the clique (persistent 
global optimum for A > A* ~ 0.61) and dense but incom¬ 
plete graphs reached as we move below A*. Because the 
clique is optimal for A > A*, anything that we measure on 
this global optimum stays constant as a function of A un¬ 
til A < A*, for which our wandering over the Pareto front 
yields a changing global optimum as A decreases. Then, 
any measurement performed on the SOO optimum will 
vary steadily with a derivative (with respect to A) differ¬ 
ent from 0. This discontinuity in the derivative indicates 
that a second order phase transition takes place (Fig. 
ic)). 

A sharp edge in the front (Fig. Bb)) indicates yet 
another second order phase transition similar to that de¬ 
scribed in Fig. Bb): The slope of the front is well defined 
as we tend towards the sharp edge from the left (yielding 
a slope (T+ such that A"*" = —(T+/(l — d+) ~ 0.3), and as 
we tend to the sharp edge from the right (now with d~ 
such that A“ = 0.01). For any A € (A“,A+), SOOs are 
well defined through equation but there are not any 
points of the front with a slope —A/(l —A) = dG ((T+, d~) 
where to locate the optimum. Instead, the same one net¬ 
work laying precisely at the sharp edge is consistently 
optimal for this range of A. Anything that we measure 
about this optimum will remain as a function of A in 
(A“, A^), but samplings of the front run smoothly below 
A“ and above A+, with well defined derivatives for any 
order parameter as a function of A. Hence, two disconti¬ 
nuities are evident in this derivative (Fig. Qc), inset). 

Qualitatively, browsing the front through A is a con¬ 
tinuous transition from the clique (which is a global op¬ 
tima for a wide range of A), through the star graph 
(also a persistent optimum for a continuum of A), to the 
MST. The weighted density of edges (ti = p®) penalizes 
large connections first, which are dropped as we leave 
the clique. But enough of them survive among Pareto 
optimal graphs so that the star can be reached continu¬ 
ously, without needing a drastic rewiring that would leave 
an imprint in the order parameters. Note that these few 
long edges survive because they enable a low average path 
length (^2 = (1)*), which is still measured as the number 
of hops between nodes. Finally, to rearrange the star into 
the MST, the surviving long edges are replaced by wind¬ 
ing branches that extend visiting many nodes on their 


a 



Order parameter (A) 


Average path length 


FIG. 5: Partly geometrical problem on a circle, (a) Again, the front follows the archetype of the topological problem 
with two roughly perpendicular stretches that trade off between the clique (top-left), star networks, and the MST. Around the 
clique it is observed the same phase transition as before. In the cavity it is solved a complex rearrangement. Networks that 
drop their larger connections first (A) must morph into a star (B), which requires some of these far-reaching edges. Therefore, 
some Pareto optimal networks are produced that never get to be SOO optima (C), yielding a first order transition, (c) Order 
parameters as a function of A reveal the hrst (A“ ~ 0.34, magnihed in the inset) and second (A* ~ 0.59) order transitions, (d) 
The landscape potential unveils the mechanisms for local equilibrium and hysteresis associated to first order transitions. At 
low levels of the control parameter (A ~ 0.1) only one minimum exists in the global energy fl. A pocket becomes locally stable 
for A ~ 0.2. This grows for larger A, until it becomes the global extreme of the energetic landscape (A ~ 0.4). Optimizing our 
networks through numerical algorithms can get us stuck in local minimums, so to transit from one potential well to the other 
we need to increase our control parameter until one of the wells get destabilized (A ~ 0.57). Repeating the operation with 
decreasing A can get us stuck in the other well, thus engaging in a hysteresis loop. Some Pareto optimal networks inside the 
cavity are reached at these metastable states. 


way. Alternative strategies, like hybrid MSTs that incor¬ 
porate non-essential shortcuts between far-apart nodes, 
fall off the Pareto front (note the sub-optimal graph high¬ 
lighted in Fig. I^a)). 

It is very useful to consider the potential landscape in¬ 
troduced by equation ([^ to stretch our intuition. For 
a fixed value of A we compute the energy 17(7, A) for 
every Pareto optimal network. The result is a lower en¬ 
ergy boundary (Fig. Qd)). Not Pareto-optimal networks 


must present yet higher energies. The SOO solution be¬ 
comes now very intuitive since the global optimum lays 
at the minimum of 17(7, which has got a vivid graphic 
representation. But this potential landscape changes as a 
function of A, and consequently its minimum. The sharp 
edge of the Pareto front is inherited by the potential land¬ 
scape, which presents a persistent minimum for a range 
of the control parameter. 
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2. Nodes spaced over a circle 

The front of this problem (Fig. [^a)-(b)) again fol¬ 
lows the archetype: i) a stepped stretch to the left that 
trades-off between the clique and dense (yet incomplete) 
graphs and ii) a flatter stretch that encompasses graphs 
with roughly the same edge density but a wide variation 
along the average path length dimension. This second 
stretch extends to large values of ti: a region populated 
by minimal, circlelike networks with little long-range con¬ 
nections. It seems a convex stretch that ends smoothly, 
so nothing remarkable happens there. Again, the stepped 
stretch of the front is convex (Fig. [^a), inset) and ends 
abruptly revealing a second order phase transition. It is 
similar to the one encountered before by the clique and 
happens at a similar value A ~ 0.59 (Fig. [11(c)). 

The notable feature that this front introduces is a 
concavity at the junction between the two archetypal 
stretches of the front (Fig. j^b)). The cavity lays at the 
confluence between the three relevant network topolo¬ 
gies: incomplete cliques, the star, and encircled nets (the 
MST of the problem). As in the previous case, longer 
edges are dropped first. The scattered nodes managed to 
retain enough long-range links in that example, enabling 
a continuous transition through the star. This is not 
possible now, and the symmetry of the circle might be 
crucial therefore. Earlier, the distribution of lengths were 
varied, while now all long-range edges are the same: the 
moment one is dropped, the others follow. As we leave 
the clique, we converge quickly to encircled graphs with 
little long range connections; and these lay inside a cavity 
of the front (Fig. Qb)). To reconstruct a star (which re¬ 
mains Pareto optimal due to its low average path length) 
a drastic rewiring in unavoidable. This prompts a first 
order phase transition (A*^ ~ 0.34) whose imprint is, in¬ 
deed, that cavity. That transition is reflected in any order 
parameter 9 that we plot as a function of A (Fig. j^c), 
inset). 

We can resort again to a potential landscape to vi¬ 
sualize this transition. Plotting 14(7, A) for all Pareto 
optimal networks we obtain the lower energy boundaries 
portrayed in Fig. [^d). This landscape changes as A 
varies, producing two potential wells associated to local 
minimums. At A = A°, both minimums present the same 
energy, thus both phases coexist. Moving away from the 
transition point, one of the wells is unstabilized. Note 
that moving A back and forth could get us temporarily 
trapped in metastable states (the most energetic local 
minimum) and hysteresis loops would be observed. 

C. Fully geometric problem 

Introducing geometry in both target functions has the 
effect of smoothing the Pareto front, removing first order 
phase transitions. Some relevant second order transitions 
disappear. Others persist, but only at the extremes of the 
front. The picture becomes closer to a soft trading-off 


between the clique and the MST. 

1. Nodes scattered over a plane 




FIG. 6: Fully geometrical problem for nodes scattered 
over a plane. (a) The front has no accidents. It is 
completely convex and spans all possible slopes so that each 
A £ (0,1) poses an SOO with a different solution. As we roll 
over the front, the clique gently leads to less connected net¬ 
works, towards the MST. (b) The absence of phase transitions 
renders smooth plots of any order parameters. 

This problem presents a quite uninteresting front (Fig. 
j^a)) without phase transitions. The front spans all pos¬ 
sible slopes d G (—00,0), so that SOOs with different 
solutions can be posed for each A G (0,1). Any order pa¬ 
rameter renders a continuous plot (Fig. j^b)) even when 
zooming in to tiny details (inset). The first derivative of 
the order parameters also behaves properly. 

All phase transitions from previous cases have van¬ 
ished. Besides the geometric disposition of the nodes 
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(which has some obvious influence over what transitions 
are present), it is notable that choices of optimization tar¬ 
gets exist for which previously existing transitions disap¬ 
pear. The fact that, given a same set of networks 7 G T, a 
choice of targets erases previous transitions implies that 
a gentle evolution between radically different topologies 
(clique and MST) can happen, despite the drastic mod¬ 
ifications that we might envision necessary a priori, and 
despite the phase transitions that do take place on the 
same graphs for other choices of targets. This stresses 
the role of target functions to frame phase transitions 
properly. 


2. Nodes spaced over a circle 

Again, introducing geometry in both targets smooths 
the Pareto front (Fig. [^a)). The first order transition 
found in the circle before has disappeared. There is no 
cavity now and the tradeoff between clique and circle 
happens gradually as we roll over the front. That previ¬ 
ous transition took place because the gradual drift from 
clique to circle was interrupted by the presence of the 
Pareto optimal star, which kept the average path length 
low because it was measured as the number of hops be¬ 
tween nodes. But now geometry also enters through 
ti = {ly ( 7 ) and using only two links to get from one 
node to another is still costly if these are far reaching 
connections. It is more economic now to circle around 
even if that implies visiting many more nodes. Thus the 
star is retracted from the Pareto front and the transition 
from the clique to the MST proceeds smoothly. 

The bottom stretch of the front seems convex but 
abruptly terminated, suggesting a second order transi¬ 
tion at Aj ~ 0.28 that trades off between the MST (an 
almost complete circle, which is persistently optimal for 
A < Aj) and other, more connected graphs. The char¬ 
acteristic plot of a second order transition is noted in 
any order parameter (Fig. Hb)). At the other end of 
the front we find the usual transition associated to the 
clique, which did not disappear but has been moved to 
A 2 — 0.98. The same characteristic order parameter plot 
can be appreciated (Fig. Bb), inset). For A G (A^jA^) 
any order parameter is a smooth function of A. 


IV. CONCLUSIONS 

In this paper we solve three MOOs defined on complex 
networks. These problems allow us to explore interesting 
aspects of Pareto optimality. Following recent contribu¬ 
tions in biology [HI 02 ] , our work is an exploration of a 
morphospace. A first approach to such spaces is to list 
all possible morphologies for a system and locate them 
quantitatively with respect to some relevant aspects - 
here, complex networks are characterized in an average- 
path-length vs. edge-density two-dimensional space. We 
propose that a natural selection process based on Pareto 




FIG. 7: Fully geometrical problem on a circle, (a) The 

front presents a smooth transition between the clique and the 
open circle, with no relevant feats except in the extremes of 
the front. These end up abruptly, as in second order phase 
transitions, (b) Plotting any order parameters reveals these 
phase transitions at A) ~ 0.28 and A^ — 0.98 (inset). 


optimality shall constrain further such a morphospace 
and we study the effects of these conflictive restrictions. 
In doing so, we follow recent works [43ll46] that illustrate 
how Pareto optimality reduces the effective dimensional¬ 
ity of certain complex systems. 

Alternatively, a recent theoretical framework brings to¬ 
gether statistical mechanics and MOO [35], thus enrich¬ 
ing the analysis of Pareto optimal designs. This frame¬ 
work reveals universal features of Pareto optimal systems 
that correspond to phase transitions or critical phenom¬ 
ena. 

First order phase transitions indicate that a system 
must undergo important structural changes despite little 
variation of some control parameter. In thermodynam- 
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ics, this implies great investments of energy in reshaping 
matter, e.g., as it transits from solid to fluid. Similar 
demands might be requested of complex Pareto optimal 
systems, specially if the parameters controlling the con¬ 
ditions for optimality may change over time. This un¬ 
derscores the importance of gathering knowledge about 
the Pareto front before implementing solutions to any 
optimization problem. 

Second order phase transitions can also be very infor¬ 
mative about the nature of a system. They indicate that 
some solution is stable for a large range of the control pa¬ 
rameters, thus making it more likely if evolution or design 
has taken place under many different scenarios. On the 
other hand, if such stable solutions would not show up 
in an evolutionary setup, we would have strong evidence 
that a large set of possible circumstances do not occur 
naturally. 

We have found a variety of first and second order tran¬ 
sitions. These depend very much on the precise mathe¬ 
matical expression of the optimization targets (see, for 
example, how transitions disappear as we change the 
measure of average path length). Hence, looking at the 
problem from an alternative perspective, the presence 
or absence of expected phase transitions in real systems 
could be informative about the nature of the optimiza¬ 
tion pressures that these systems might be subjected to. 

Finally, the fully topological problem presents a quite 
singular case: a first order phase transition takes place 
between the star graph (for A < A^) and the clique 
(A > A'^), while all other Pareto optimal networks are 
also global optima at A = A*^. This would suggest that 
cliques and star graphs should happen overwhelmingly 
more often than any other topology when geometry is not 
relevant, unless every optimal network had evolved under 
the unlikely condition A = A'^. This is notably at odds 
with the reality and a possible solution (with connections 
to critical phenomena) will be explored elsewhere |40j . 
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Appendix A: Analytic and numeric approaches for 
MOO solving 

We relied on genetic algorithms to locate the different 
Pareto fronts. The fully topological case can be solved 
analytically, but the same genetic algorithm was used to 
check for good convergence showing very good results. In 
this appendix we explain in detail the genetic algorithm. 
Following that explanation, a few disclaimers are in or¬ 
der about the numerical nature of the solutions found - 
i.e. about the fact that convergence to the Pareto front 


cannot be guaranteed and the smoothing necessary to 
render a continuous approximation to the front. 


1. A multiobjective genetic algorithm 

MOO relays on the concept of Pareto dominance. 
Given two networks 7 i, 7 j G F, both mapped into 
through ti( 7 i/j) and ^ 2 ( 71 /^), we say that 7 ^ dominates 
7 _; (and note it 7 ^ ^ 7 ^) if 7 ^ is not worst than 7 j re¬ 
garding any target and it is better than 7 j with respect 
to at least one target. We can visualize this: Since we 
deal with minimizations in the ti — plane, network 
7 i has got a set of axes associated with their origin at 
(ti( 7 z),^ 2 ( 7 *)) and every network 7 ^ laying on the first 
quadrant of these axes is Pareto dominated by 7 ^. 

Following the literature on multiobjective genetic al¬ 
gorithms [17IH5] we computed a dominance score: We 
took the set Dj = { 7 i| 7 i ^ 7 ^} of solutions from within a 
given population (an arbitrary subset of F) that dominate 
7 j. The size of this set (dj = dubbed the dom¬ 

inance score), indicates how fit 7 ^- is in terms of Pareto 
optimality. We proceeded then to minimize this score. 
We departed from an initial population of JVp networks 
(either random or designed, see below), selected Np/2 of 
the population based on the dominance score, chose ran¬ 
dom pairs among the selected networks to produce Np/2 
new networks, and applied random mutations to all but a 
subset of elite networks. We iterated this scheme a fixed 
number of generations. 

Mutations consisted in random appending or deleting 
edges or totally swapping the connections of two arbi¬ 
trary nodes. For crossover, from each of the two mating 
nets we assigned each node and its connections randomly 
to each of the offspring graphs checking that the same 
node and connection was not assigned twice to the same 
child. After crossover or edge deletion we checked that all 
networks remained connected all the time. We completed 
one missing link whenever connectedness failed and then 
checked for connectedness again. 

All we care about for the current research is good 
convergence towards the front. This justifies our using 
of crossover: this is a very good evolutionary operator, 
though unrealistic if we wanted to study some features 
of nature. For example, such an operator would not be 
adequate to study species that do not reproduce sexu¬ 
ally. Studying Pareto optimality under constrained con¬ 
ditions - e.g. without crossover - also renders a set of 
non-dominated solutions. These might converge to the 
Pareto front or not, and they might be subjected to ge¬ 
ometrical constraints in ti — t 2 that are similar to those 
studied in |35j and in this paper for Pareto optimal net¬ 
works. Such constrained evolutionary schemes pose inter¬ 
esting research questions, but here we are concerned with 
Pareto optimal solutions. This justifies the crossover and 
a clever initialization of the algorithm. It might be diffi¬ 
cult to converge towards some solutions that are highly 
non-trivial - e.g. the MST. We know, though, that this 
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solution belongs to the front of all physically grounded 
problems. If an algorithm would fail in hnding this so¬ 
lution, this could hinder convergence towards an inter¬ 
esting (though challenging) region of the Pareto front. 
Once again, because we are concerned with Pareto op¬ 
timal solutions and would like to attain the closest con¬ 
vergence possible, it is also fully justihed to seed the ini¬ 
tial population with a few designed solutions. We did 
so by introducing since the very beginning cliques, star- 
graphs, MSTs, and circle networks (these two are equiv¬ 
alent for the circle) with very slight mutations. We pro¬ 
duced 7Vp/4 of each such major topologies at the begin¬ 
ning. The crossover and mutation operators ensure fast 
exploration of hybrid topologies. 

As noted above, we know that some of these networks 
are Pareto optimal: the clique is always so, and the MST 
is Pareto optimal in all physically grounded problems. 
However, and since we started with little variations upon 
these graphs, the algorithm did not always reach these 
solutions - but it surely explored the region nearby. 
There might be other interesting regions of the front 
that might not have been fully explored and that are 
impossible to seed without foreknowledge. Although we 
are concerned with Pareto optimal solutions alone, our 
methods are numerical at the end and convergence of 
multiobjective genetic algorithms to the Pareto front 
cannot be guaranteed. We decided to report on the 
results of the simulations with as little reinterpretation 
and further speculation as possible. Notwithstanding, 
the overall details of the Pareto front seem to be recov¬ 
ered and the theory posed in is properly illustrated. 

As for the implementation of the algorithm, we used 
a population of Np = 3000 connected networks with 
N = bO nodes - with the initial population seeded as 
indicated above. The population was evolved during 
T — 10000 generations in every case. Mutation hap¬ 
pened with a probability = 0.001 of appending an ex¬ 
tra link to each network, the same probability of deleting 
an existing edge, and the same probability of swapping 
the ends of each existing connection. The top Ng. = 50 
networks of the population where considered elite and 
were spared any mutation. As the algorithm proceeds, 
many networks reach a dominance score of 0 even if they 
are not Pareto optimal. Unluckily, this score is the best 
indicator of Pareto optimality available (not only in the 
current implementation, but generally). This results in 
elite members of the population not being objectively 
better than non-elite members - in terms of Pareto opti¬ 
mality. Because the algorithm sorted the population con¬ 
sistently from one generation to the next, what members 
of the population are considered elite is largely a matter 
of antiquity: early members that reach low dominance 
score and are not overthrown are likely to be preserved 
during the whole simulation. We repeated 4 times the 
simulations with scattered nodes to check that the rele¬ 
vant features obtained were not artifacts of some lucky 
distributions of the nodes. 

2. Smoothing of the front and order parameters 

The only speculation that we allowed ourselves is in 
choosing a relevant scale for analysis, which led to a 


smoothing of the Pareto front. As noted in Sec. Ill 


we deal with a discrete set of networks whose front is 
necessarily discrete as well. Accordingly, every shift in 
global optima is a hrst order phase transition at some 
scale and global optima remain so for a continuous range 
of A, as in second order phase transitions. This does 
not further our understanding of the problem as much as 
a coarse-grained analysis that renders noteworthy phase 
transitions. Since the genetic algorithm only produces 
a finite set of (ideally) Pareto optimal solutions, we ap¬ 
plied a Bezier smoothing to their plot on the ti — t 2 plane. 
We took care that the smoothing did not introduce alien 
concavities. Because Bezier curves cannot present sharp 
edges (thus ruling out second order phase transitions), 
when a sharp edge seemed the best description of the 
front (Sec. IIIB partly geometrical problem on nodes 
scattered over a plane), we decided to split the front in its 
two salient branches and apply two independent smooth¬ 
ing processes that allowed us to recover the transition in 
great detail. 


To locate global optima, we calculated H(a;n, A) for 
the optimal solutions produced by the genetic algorithm, 
and for a large sample of points from the Bezier curves 
introduced in the previous paragraph. We registered 
the global optimal for different values of A. One of 
the problems pointed at earlier is that, because of the 
discreteness of the front, global optima are so for several 
values of A. This would cause that the plots of order 
parameters look tiered. For a better illustration of 
the results, whenever order parameters are plotted we 
indicate only the first and last values of A for which 
each global optima are indeed optima (black crosses 
in all order parameter plots). The smoothing allows a 
finner grained sampling so that this is not an issue: the 
corresponding order parameters (red curves in all order 
parameter plots) look continuous always. 


Following |35j . anything well behaved that we measure 
upon global optima are accepted as order parameters. By 
well behaved we imply that order parameters should not 
introduce alien divergences into the problem, and that 
solutions laying at different points over the front should 
score differently in this parameter. This way we ensure 
that any feature stemming from the optimization prob¬ 
lem does not go unreported and that we do not introduce 
phase-transition-like behaviors that originated, e.g., on 
some function diverging to infinity for reasons of its own. 
Taking these guidelines into account, the target functions 
themselves are always good order parameters. We use 


these {0 — ti in Sec. Ill C 

, or trivial transformations of 

them (0 = 1 — ^2 in Sec. 

IIIB 

). More drastic transfer- 


mations such as 1/(2 — ti) would be banned: note that 
this function diverges for ti = 2 even if this is a perfectly 
regular point of the front for all problems. 
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